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ABSTRACT 

The galactic central black hole Sgr A* exhibits outbursts of radiation in the near 
infrared (so-called IR flares) . One model of these events consists in a hotspot orbiting 
on the innermost stable circular orbit (ISCO) of the hole. These outbursts can be used 
as a probe of the central gravitational potential. One main scientific goal of the second 
generation VLTI instrument GRAVITY is to observe these flares astrometrically. Here, 
the astrometric precision of GRAVITY is investigated in imaging mode, which consists 
in analysing the image computed from the intcrfcrometric data. The capability of the 
instrument to put in light the motion of a hotspot orbiting on the ISCO of our central 
black hole is then discussed. We find that GRAVITY'S astrometric precision for a 
single star in imaging mode is smaller than the Schwarzschild radius of Sgr A*. The 
instrument can also demonstrate that a body orbiting on the last stable orbit of the 
black hole is indeed moving. It yields a typical size of the orbit, if the source is as 
bright as mx = 14. These results show that GRAVITY allows one to study the close 
environment of Sgr A*. Having access to the ISCO of the central massive black hole 
probably allows constraining general relativity in its strong regime. Moreover, if the 
hotspot model is appropriate, the black hole spin can be constrained. 

Key words: Instrumentation: interferometers - Astrometry - Galaxy: centre - Black 
hole physics 



1 INTRODUCTION: TESTING GR IN ITS 
STRONG REGIME WITH GRAVITY 

GRAVITY is a second generation VLTI instrument the 
main goal of which is to test general relativity (GR) at the 
very centre of the Galaxy (e.g. Eisenhauer et aL||2008 



Table 1. Main characteristics of the GRAVITY instrument. 



also Table [T]). Observing the motions of sources orbiting in 
the vicinity of the dark mass Sgr A* with great accuracy 
will allow probing the gravitational potential of the central 
massive black hole. 

The star observed closest to Sgr A* so far, S2, has been 
monitored for more than ten years, and its orbit is known 



2002 Ghez et al. 2008 



with high precision (Schodel et al 
Gillessen et al.|2009 |. This allowed to constrain the nature 
of the central attractor and to advocate with a high level 
of confidence that it is a black hole. However, the 15.8-year 
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Maximum baseline length 

Number of telescopes 

Aperture of each telescope 

Wavelengths used 

Angular resolution 

Size of the totaj^j field of view 

Size of the scienccj^] field of view 



143 m 

4 (all UTs) 

8.2 m 

1.9 - 2.5 fim 

4 mas 
2" 

71 mas 



a Containing the science target and a phase reference star 
^ Containing only the science target 

period of S2 does not allow to test relativistic effects on 
the orbit within a few years of observation. For instance, 
the advance of the peribothron will probably not be de- 
tectable before the next pericentre passage (so after nearly 
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20 years of cumulated observations). A powerful test of GR 
therefore needs to have access to stars closer in, at least ten 
times closer to Sgr A* to be able to measure orbits with 
periods of the order of a year ( |Rubilar &; Eckart|2001 Pau 



netic reconnection orbiting on the innermost stable circular 



mard et al. 2008 1 for which the relativistic effects will be 
greater and detectable over a shorter period of time. Such 
stars, fainter than S2 and more concentrated toward Sgr A* 
than the currently known S-stars, will be more confused: the 
highest possible resolution is then crucial to observe their 
orbital traces. Will (20081 shows that a 10 /xas astrometric 



precision (the goal of GRAVITY) is sufficient over one year 
to measure peribothron advance, frame dragging and the 
quadrupole momentum, provided stars with short enough 
periods and high enough eccentricities are available. Such a 
measurement would allow testing the 'no-hair' theorem of 
GR. 

GRAVITY aims at detecting these innermost sources 
in its field of view of 71 mas FWHM diameter. The analysis 
of the galactic centre stellar population ( Genzel et al.|2 003b 



Paumard et al. 2006), extrapolated to the central 71 mas, 



allows one to estimate the number of stars that GRAVITY 



will be able to detect. Fig. 7 oflGenzel et al. (2003b I gives a 



number density of around 25 stars with 10 < 771k < 17 per 
square arcsecond with angular distance to the centre less 
than 35 mas (i.e. inside the field of view). There are thus 25 x 
tt x 0.035 2 w 0.1 star in GRAVITY'S field of view. Thus, less 
than one star in this range of magnitude is expected to be 
continuously detectable (as the density function is increasing 
with decreasing radius, more fainter stars will be detectable, 
but will demand a longer integration). However, stars from 
further out will fly through the central 71 mas from time 
to time due to their large enough eccentricities, such as the 
three mK = 15 S-stars: S2, S14 and S17 nowadays. Thus, 
considering around 3 stars of magnitude iuk ~ 15 to be 
in the field of view is slightly optimistic, but not unrealistic. 
This will be the case, provided that the period of observation 
is well chosen (in order for the closest S-stars to be present). 

Sgr A* also exhibits powerful flares in the X-ray, NIR 
and sub-mm domains (e.g. Baganoff et al.|2001 Ghez et al 



2004{|Clenet et al.|2005||Eckart et al.|2009| ). The NIR events 
are characterized by an overall time-scale of the order of 
one hour, and a quasi-periodicity of roughly 20 minutes; 
the luminosity of the source increases by a factor of 20 (see 
Hamaus et al.|[2009). The typical flux density at the max- 



imum of the flare is of the order of 8 mjy (Genzel et al. 
2003a, Eckart ct al. 2009) which translates to approximately 
mK = 15. The brightest flare observed to date reached 
thk = 13.5, the median flux of the source thus increasing 
by a factor of 27 ( |Dodds-Eden et al.|2010 1 



It must be noticed that some authors question the fact 
that flares can be understood as specific events. The light 
curve fluctuation is then interpreted as a pure red noise (see 
e.g. Do et al.|2009 1. Following this rationale, it is no longer 
justified to select the 'flare' part of the light curve: the quasi- 
periodicity is then no longer significant. This debate of the 
event nature of the flares is investigated in|Dodds-Edc n~et al.| 
(20101. 



Various models are investigated to explain these flares: 
adiabatic expansion of a synchrotron-emitting blob of 
plasma ( Yusef-Zadeh et aL||2006 I , heating of electrons in a 
j et (|Markoff et al.|2001[ ), Rossby wave instability in the disc 
( Tagger & Mclia 2006 ) , or a clump of matter heated by mag- 



orbit (ISCO) of the black hole ( 


Yuan et al.||2004 Paumard 


et al.||2008 Hamaus et al.||2009 


Eckart et al.||2009l. To be 



tested, this last model requires an astrometric precision at 
least of the order of the angular radius of the ISCO, which 
is a few times the Schwarzschild radius of the black hole, i.e. 
a few times 10 fias. 

GRAVITY is a beam combiner that will interfere the 
light of the four VLT UTs. It is designed to cope with such a 
challenging measurement as the follow-up of an object orbit- 
ing on the ISCO of Sgr A*. A crucial requirement to do this 
is to be able to correct properly and in real time the atmo- 
spheric turbulence. Turbulence makes the phase jitter and it 
is necessary to define a phase reference and to maintain the 
fringes centered on this reference. This will be made possi- 
ble by observing an unresolved (with zero intrinsic phase) 
reference source at the same time as the science target. The 
phase of the reference being null, it can be used to define the 
zero of the optical path delay. A dedicated fringe tracking 
system is used in GRAVITY to keep the fringe around the 
zero level (for details, see |Gillessen et al.|2010 1 . 

It appears clear that astrometric precision is a key el- 
ement to the realization of the GR-related science cases of 
GRAVITY. As the instrument will be phase referenced, will 
use the four VLT UTs and will have five spectral channels 
in low-resolution mode, the u-v plane coverage will be large 
enough and the precision on the phase sufficient to con- 
sider using GRAVITY in imaging mode, i.e. to reconstruct 
images from the interferometric observables. Here, it is as- 
sumed that the scientific target is made of simple continuum 
sources in order to be able to combine the various spectral 
channels. 

Determining the astrometric precision in imaging mode 
for several fixed stars present in the instrument's field of 
view is the first aim of this article. A few fixed stars (1 to 
3) are placed in a square box of size 100 mas. The interfero- 
metric data acquisition is simulated, realistic noise is added, 
the image is then computed from the noisy visibilities and 
phases. The CLEAN algorithm is used to deconvolve the im- 
age. A Monte Carlo analysis based on this procedure will 
lead to the determination of the astrometric precision. This 
part of the study is an extension of an earlier work by |Pau-| 



mard et al. (20081 



The second objective is analysing the case of a variable 
source present in the field. The typical light curve of an IR 
flare is modeled and the source is placed on a circular orbit 
of radius 30 fias (ISCO of a Schwarzschild black hole of the 
mass of Sgr A*). It will be determined wether its movement 
is detectable, and to what extent. 



2 STANDARD TWO-FIELD 

INTERFEROMETRIC ASTROMETRY VS 
IMAGING MODE ASTROMETRY 

2.1 Standard astrometric mode 

The usual method in interferometry to determine astromet- 
ric information is to record the interferograms of a refer- 
ence source and an object of interest (narrow-angle, dual- 
feed astrometry). If the reference and the object are sup- 
posed to be isolated (no other sources) and point sources, 
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the phase shift between the two interferograms is simply 
equal to: A(f> = A<^ p i a ton + ^f-B ' ^ being the baseline 
vector, AS the angular separation on sky and A<£ p i a ton the 
difference between the atmospheric piston phases affecting 
the two sources (which is not necessarily close to zero unless 
a fringe tracker is used). This is illustrated in Fig. [lj where 
the piston phase is not taken into account. Integrating for a 
time longer than the coherence time of the turbulence, the 
atmospheric piston will get averaged out and it will be possi- 
ble to deduce the separation between the sources. However, 
this method is strongly dependent on the two conditions of 
point-like and isolated sources. If one of the sources is ex- 
tended, some confusion may arise in the interpretation of the 
astrometric measurement as the phase information may be 
contaminated by the source structure, which has a non-zero 
intrinsic phase. Moreover, a fixed source could give rise to a 
changing astrometric measurement if its surface brightness 
distribution evolves with time. Even worse, the motion of 
the source components could be undetectable if the object's 
geometry is perfectly symmetric so that its flux barycentre 
is still (such as for a perfectly balanced binary), because 
the measured phase traces the barycentre position only. It 
should be noted that there is no definitive evidence (only 
reasonable assumptions) as to whether the environment of 
Sgr A* is complex or not as this will be discovered by GRAV- 
ITY in the near infrared with possibly first hints by PRIMA 
( Delplancke|2008[ ) and ASTRA ( |Pott et al.|2008 l. 

This standard method can still be used in more complex 
contexts by resorting to model-fitting. The complex visibil- 
ity is computed assuming a given model described by some 
parameters (e.g. the astrometric positions of the sources and 
their relative intensities). Fitting these modeled visibilities 
to the observed ones allows one to determine the param- 
eters. However, if the field is too complex, the number of 
free parameters will make it extremely difficult to use this 
procedure. 



2.2 Imaging mode 

Here, the 'imaging' mode of the instrument will be investi- 
gated. 

The aim is to reconstruct an image from the complex 
visibilities, using the Zernike-van Cittert theorem. It is thus 
necessary to have access to the source phase with high pre- 
cision. This is made possible by using a phase reference, 
defining the zero of optical path delay. This procedure gives 
rise to specific errors on the phase, that are described in 
detail in Appendix |A3| (see in particular equation ( |A10[ )). 

The performance of the technique relies also on the abil- 
ity to reconstruct a deconvolved image from the rough data 
obtained by applying the Zernike-van Cittert theorem. If 
the field contains only point sources, the CLEAN algorithm 
(see Hogbom et al.|1974| | will be efficient for the deconvolu- 
tiorQ The astrometric information is then directly available 
on the deconvolved image. 



2.3 Imaging mode or fitting? 

The first aim of this work is to investigate the astromet- 
ric precision of such an instrument as GRAVITY in imag- 
ing mode and in standard fitting mode. We are interesting 
in knowing whether the pure imaging mode can be of any 
scientific interest, or whether a fitting procedure is always 
necessary in order to obtain high enough precision. 

Should the imaging mode be able to give precise enough 
results, it could appear as a very interesting alternative to 
fitting. Indeed, if the astrophysical source is too complex to 
allow for a few-parameter model to be determined, astro- 
metric parameters could still be measurable directly on the 
image. Moreover, an image better catches most of the ob- 
ject's high-frequency details ( Haubois et al.|2009 1. However, 
testing the imaging mode with a complex field of view is 
beyond the scope of this work as it would need a far more 



sophisticated algorithm than CLEAN, such as WISARD (Mei- 
mon et al"]|2009| ) or mira ( Thiebaut||2008 1 . A discussion of 



the use of these two algorithms, in comparison to standard 
parametric model-fitting, can be found in |Le B csncrai s~et al.| 
(20081. 



1 Here, a Yorick implementation of CLEAN is used, which 
can be found at the URL http://www.lesia.obspm.fr/perso/ 
thibaut-paumard/yorick/CLEAN . i 



Should a fitting procedure be necessary in order to ob- 
tain the desired precision, imaging mode could still be very 
interesting to determine the broad lines of the field of view. 
This information could then be used as a first guess for the 
fitting. 



2.4 Limiting precision 

Whatever the method used, it is possible to compute the 
limiting precision that one can hope to obtain, which cor- 
responds to the theoretical precision of the instrument. The 
phase error associated with an integration time of 100 s on 
a niK = 12 source has been determined by simulations: 
CT0 = 0.0085 rad. Then, for 7V= 30 u-v points (equal to 
the number of wavelenghts - 5 - multiplied by the num- 
ber of baselines - 6 -), a mean baseline of _B=100 m, and 
a mean wavelength of A — 2.2 x 10 -6 m, the astromet- 
ric error is of: <7theo = ~ XA2 uas. This estimate 
assumes that all the u-v points are aligned in one given di- 
rection. A rough estimation of the error in a given direction 
is: <7 y thco ~ v2ftheo = 1-58 pas. 

Thus in any case, the precision will stay above approx- 
imately 1.6 /ias. 



3 SIMULATING THE DATA 

The astrophysical object of interest, a star field, is described 
by a set of coordinates in mas (2;f odel , yf ode ), randomly 
placed in a 100 mas square box, and by intensities h. The 
(u, v) plane coverage is computed and the modulus and 
phase of visibilities are derived by means of the Zernike- 
Van Cittert theorem. All the different sources of noise intro- 
duced in Appendix [X] that affect the visibility modulus and 
phase are then added, as a function of the source magni- 
tude. The noisy visibilities are then transformed back to an 
image. This is performed by computing a numerical Fourier 
transform: 
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Star 1 Star 2 




opd 

Figure 1. Standard two-field interferometric astrometry illustrated in one dimension. The baseline is B, the angular separation on sky is 
AS, and the optical path difference (opd) is related to the phase difference through: A<f> = ^ Aopd, where A is the radiation wavelength. 



fc=l 



(1) 



where the bold quantities are 2D arrays (I being the image), 
Ve lv is the complex visibility, and iV uv is the number of pairs 
of spatial frequencies (u,v). 

Because of the use of single-mode fibers, the transmis- 
sion of the 71 mas FWHM field of view is gaussian: stars 
far from the centre will be damped. Images are deconvolved 
using the CLEAN algorithm, which concentrates the flux of 
each star in a set of a few pixels of coordinates (xk,Vk)- 
The deconvolved position (x decom , y decom ) associated with 
each set of pixels is defined by computing the barycentre 
of the pixels weighted by the intensity. The typical num- 
ber of illuminated pixels varies from 2 to 12 (depending on 
the number of sources and integration time) in the different 
cases considered in this work. Thus: 



x k Ik 

k 



(2) 



v 



where Ik is the intensity of pixel number k and 7 to t denotes 
the total intensity of the pixel set. 

Fig. [2] shows an example of this procedure. 

It is then easy, in the frame of a Monte Carlo analysis 
based on this concept, to measure an astrometric precision 
which is defined as the root mean square (RMS) of the dif- 
ference between the two sets of coordinates: 



nun/ model dcconv\ 

ff x = RMS(xj — Xi ) 



(3) 



a y = RMS(l/ t m 



Vi 



All these simulated data are computed using the follow- 
ing parameters. 

We consider one night of observation chosen on 2010 
April 15. Choosing the 15 th day of the month is purely ar- 
bitrary. However, the choice of the month has an impact 
on the total time during which the Galactic centre (GC) is 
observable. April is quite a typical choice as the minimum 
(non-zero) observable length is in February and the maxi- 
mum in June. During the chosen night, the GC is observable 



for a little more than 5 hours (as a comparison, it would be 
observable for 9 hours in June). 

One elementary integration block lasts 100 s. Two inte- 
gration times are used: 5 h (i.e., the whole night), or a sin- 
gle 100 s block. When the observation is performed during 
5 h, the observing procedure is chosen as follows. Integra- 
tion on the scientific target is performed during 1 h, divided 
into elementary integration blocks of 100 s, then 30 min are 
spent for calibration, and so on during the night. The Earth 
rotation during one elementary block is supposed to be neg- 
ligible because the trace in the u-v plane during the 100s 
of integration is perfectly linear and can be replaced by its 
average. Six baselines are used (with the four UTs) and five 
wavelengths (2, 2.1, 2.2, 2.3 and 2.4 fim). 

Fig. [3] shows the u-v plane coverage obtained after 5 h 
of observation. The u-v plane has been rotated so that the x 
and y axes in the image plane are aligned with the big and 
short axes of symmetry of the PSF respectively. 

The sampling (number of pixels per mas) is chosen in 
order to keep computing time reasonable, and to obtain the 
best possible precision. When the whole night is used to 
integrate, a 3 pixels per mas sampling is chosen. Finer sam- 
pling would demand too long a computing time in the Monte 
Carlo approach. When the integration time is 100 s, a 10 
pixels per mas sampling is reasonable in terms of comput- 
ing time, and can allow getting a finer precision than with 3 
pixels per mas. However, this finer sampling is only of inter- 
est when a single star is in the field. If more than one star 
is present, the PSF of the sources overlap and contaminate 
each other, leading to the apparition of spurious secondary 
peaks. In these cases, the 10 pixels per mas sampling is no 
longer interesting. Thus, a 3 pixels per mas sampling is used 
in order to gain computing time. 

It must be noticed that these samplings fulfill Shannon's 



criterium (Shannon 1949), as the maximal frequency of the 
signal is 



1.5 mas. 



4 GRAVITY ASTROMETRIC PRECISION IN 
RECONSTRUCTED IMAGES 

In this section, the astrometric precision in reconstructed 
images is computed when one to three stars are present in 
the field. Magnitudes are varied in the range 12 ^ ttik ^ 17. 
The integration time will be equal to 5 h (a configuration 
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Simulated observation oi 3 stars, itik=15, whole night integration 



Associated CLEANed image 
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Figure 2. Top left panel: image of a simulated observation of three stars with integrated magnitude ttt-k = 15 made over one complete 
night. Top right panel: illuminated pixels obtained after deconvolution of the image by the CLEAN algorithm, the blue triangles show 
the positions of the sets of pixels. Bottom panel: zoom on one of the sets of pixels. 

Whole night integration PSF u-v plane 
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Figure 3. PSF central peak and u-v plane 



which can be used in practice to image a field of nearly 
motionless stars), or to a single 100 s block (which can be 
used to image a quickly evolving field, such as a flare). 

For each Monte Carlo simulation made of N different 
runs, run number 1 ^ i ^ N gives a value of the differences 
df and d v . between deconvolved and model coordinates: 



d, 



model deconv 
Xj Xj . 



(4) 



These quantities are in turn used to compute the RMS 
defined in equation |3l. 

The one-star case is of purely theoretical interest in or- 
der to determine the best possible precision in reconstructed 



images. Indeed, when the field of view contains only one 
point-like source, the standard two-field narrow angle astro- 
metric method would provide better results because para- 
metric reconstruction (i.e. fitting) is more efficient than non 
parametric methods (such as our imaging mode algorithm) 
as soon as the model is correct (see Le Besnerais et al.|2008 1 . 
However, we will show that in the 100 s integration case, the 
imaging mode results are very close to the fitting ones. 



Figs. |4j [5] and [6] show the resulting astrometric preci- 
sions for the two different integration times and in the two 
directions x and y. 
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Astrometric error with one star, whole night integ. Astrometric error with one star, 100s integ. 
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Figure 4. Imaging mode astrometric precision in ^tas when one star is present in the field of view as a function of the magnitude; in 
the x direction (blue) and y direction (red). The total integration time is equal to the whole night (left) or to 100 seconds (right). The 
solid horizontal lines show the level of precision without any noise added to the complex visibilities. 

Astrometric error with 2 stars, whole night integ. Astrometric error, brightest ot 2 stars (+2 mag), 100s integ. 
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Figure 5. Imaging mode astrometric precision in /xas when two stars are present in the field of view as a function of the integrated 
magnitude; in the x direction (blue) and y direction (red). The total integration time is equal to the whole night (left) or to 100 seconds 
(right). In this last case, only the brightest source is searched for, which is two magnitudes brighter. The solid horizontal lines show the 
level of precision without any noise added to the complex visibilities. 



4.1 General remarks 

The astrometric error as a function of magnitude can be 
described with two main components: a plateau up to a cut- 
off magnitude and a steep increase of the error at fainter 
magnitudes. The plateau corresponds to a regime in which 
the astrometric accuracy is dominated by numerical errors 
of the CLEAN procedure as shown by the horizontal solid 
lines in Figs. |4j [5] and [6] These solid lines show the level 
of precision obtained with noiseless data, i.e. when the only 
error is due to CLEAN (which will be labeled 'CLEAN error' 
in the following). 

The steep increase of the error at fainter magnitudes is 
due to the increasing importance of the detection noise. 

The value of the plateau-limited precision is directly 
linked to the sampling. For instance in Fig. [4j the plateau is 
at 12 fias for the 5 h case (3 pixels per mas sampling), and 
3.7/ias in the 100 s one (10 pixels per mas sampling). Would 
the sampling be set to 3 pixels per mas in the 100 s case, 
the plateau value would rise to approximately 12 fia.s, as for 
the 5 h case. 

Moreover, the better precision obtained in the y direc- 
tion is due to the shape of the PSF the axes of which are not 
equal: the y axis being shorter, the precision in this direction 
is better. 



We have tried to use one of the sources in the field of 
view as a reference in order to determine the relative astro- 
metric precisions on the other two stars. However, even with 
a purely correlated error, the simulations show that there is 
no gain in using a position reference inside the field of view: 
the reconstruction errors in the coordinates retrieved are too 
high and spoil the relative precisions. 



4.2 Skipped Monte Carlo runs 

For some simulations, when the noise level is big enough, 
the CLEAN algorithm can behave badly and generate out- 
put coordinates very different from the model ones. These 
'bad runs' are very rare and can be neglected, not to pollute 
the resulting precision. Indeed, only a very small fraction of 
such 'bad runs' can drive the precision from a few tens of 
microarcseconds to the arcsecond. 

In order to decrease the number of 'bad runs', the stan- 
dard procedure described in Section [3] is adjusted for some 
cases. 

When the integration time is equal to 5 h, 'bad runs' 
occur when more than one star are present in the field and 
when they are very close to each other. Thus in these cases, 
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Astrometric error with 3 stars, whole night integ. Astrometric error, brightest ot 3 stars (+2 mag), 100s integ. 
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Figure 6. Imaging mode astrometric precision in /jas when i/iree stars are present in the field of view as a function of the integrated 
magnitude; in the x direction (blue) and y direction (red). The total integration time is equal to the whole night (left) or to 100 seconds 
(right). In this last case, only the brightest source is searched for, which is two magnitudes brighter. The solid horizontal lines show the 
level of precision without any noise added to the complex visibilities. 
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Figure 7. Histogram of the quantities d x (blue) and d v (red) defined in equation |4jl with one hik = 17 star in the field and 100 s of 
integration. Besides the gaussian shaped distribution (left), three realizations occur at values of d x and d y between 2 and 10 mas (right), 
thus disconnected from the broad distribution. These are 'bad runs': due to a high level of noise in the image, the sharpest peak of the 
image is shifted from the true position of the source. 



the stars are forced to be further from each other than 5 
mas, the PSF central peak's characteristic dimension. 

When the integration time is equal to 100 s, there are 
two main ways to decrease the number of 'bad runs'. First, 
the model coordinates (:r modo1 , j/ modo1 ) can be restrained to 
be picked in a disc of radius 35 mas (i.e., the gaussian field 
of view of GRAVITY) instead of being picked in the whole 
100 mas square box. Second, the CLEAN algorithm itself 
can be somewhat modified in order to prevent pathological 
behaviors. The standard working of the CLEAN algorithm is 
the following. 

Let us consider an observed image and a 'CLEANed' 
image, initially identically zero. 

(i) Step 1: The maximum value of intensity in the ob- 
served image is located. 

(ii) Step 2: The same pixel in the CLEANed image is added 
a given fixed value of intensity. 

(iii) Step 3: The corresponding amount of the PSF, mul- 
tiplied by a loop-gain factor, is subtracted at the maximum 
location on the observed image. 

(iv) Back to Step 1. 

The stopping condition is that the new maximum of 



the observed image becomes less than a given fraction of its 
initial maximum: 
maxinit 



maXnew ^ 



D 



where D is a chosen parameter, the value of which is taken 
between a few and 10. 

Problems can arise when both the noise and D are 
high. Then, by subtracting many PSF on the true maxi- 
mum of the image, some spurious secondary maxima may 
become sharper in the CLEANed image, leading to a com- 
pletely wrong retrieved value of the coordinates. An easy 
way to prevent this is to verify at the end of the CLEAN 
algorithm that the CLEANed image maximum is close to the 
initial maximum of the observed image (of course this only 
insures that the brightest star is not a fake in the CLEANed 
image, but we will see in the following that only the brightest 
star can be retrieved when many sources are present with an 
integration time of 100 s). If these two maxima do not coin- 
cide, then the CLEAN algorithm is restarted with a smaller 
value of D. 

These modifications of the standard procedure defined 
in Section[3]are always used (for an integration time of 100 s) 
when more than one star are present in the field, and even 
with a single star, when tuk ^ 16.5. 
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However, despite all these precautions to prevent 'bad 
runs', some may remain: if, due to noise, the maximum of 
the observed image is shifted from the true star position, 
there is no way to find the source's position with CLEAN. 
These pathological cases will be removed for the final pre- 
cision computation. In the following sections, it is always 
specified whether 'bad runs' were removed, and the percent- 
age of these skipped runs is provided. 



4.3 One-star case 

Here, 'bad runs' only occur when the integration time is 
100 s and for m,K = 17. Fig. [7] shows the histogram of the 
differences df and d\ defined by equation Q in this case. 
It appears that, besides a peak located near zero, some rare 
events are present at quite higher values. In order not to 
pollute the result, the final precision is computed without 
using the runs that end up with df or d\ greater than 2 mas 
(the first polluting peak in Fig. [7] is indeed located a little 
higher than 2 mas). The proportion of runs skipped in this 
procedure is of 0.03 per cent. 

The CLEAN algorithm is used here in its optimal con- 
figuration with only one point source. In this specific case, 
its result is nearly independent of the shape of the PSF: the 
CLEAN error is weakly dependent on the direction and on 
the integration time, it is mainly dictated by the sampling, 
as stressed in Section [4.11 

Of course, the precision cannot be ever more improved 
even if the sampling is refined: it will always remain greater 
than the theoretical limiting precision. 

The best precision obtained is by far better than the 
pixel size. For instance in the 5 h integration case, the pixel 
size is 0.33 mas so the precision is better than the pixel size 
by a factor of 25. This is due to the fact that the barycen- 
tre of the illuminated pixels is calculated, which makes it 
possible to obtain sub-pixel precision. 

It must be emphasized that the best precision obtained 
in the 100 s case, 4 /xas, is very close to the theoretical 
limit: this is a powerful demonstration of the interest of the 
imaging mode. Moreover, this theoretical precision can be 
obtained provided the sampling is high enough. With 150 
pixels per mas, the pure imaging mode precision in the y di- 
rection reaches 1.6 fias. However, this is done at the expense 
of computing time. 

It is interesting to check whether one reaches the theo- 
retical limit error when using the pure imaging mode results, 
with a 10 pixels per mas sampling, to bootstrap a fitting pro- 
cedure. The complex visibilities are fitted to the simulated 
data, the only parameters being the x and y angular posi- 
tion of the source on the sky. The intensity is not fitted, it 
is supposed to be measured with a 10 per cent error by the 
adaptive optics of the instrument. 

The imaging mode is used to determine a first guess of 
the position of the star in the field. More precisely, the first 
guess is chosen by sampling a box centered on the imaging 
mode retrieved position, its dimension being of the order of 
the imaging mode precision at the considered magnitude. 
The iteration is stopped as soon as the % 2 corresponding 
to a position in the box is between 0.5 and 1.5. Then the 
fit is performed with Yorick's lmf it routine. This procedure 
allows to have access to a precision of 1.6 /xas when inte- 



grating 100 s on a % = 12 source. This coincides with the 
theoretical value given in Section [2. 3| 

4.4 Two-star case 

With two stars in the field and when integrating during the 
whole night, Fig. [5] shows that the astrometric precision is 
of the order of the Schwarzschild ISCO angular radius (30 
pas) in both directions and for all magnitudes. 

When the integration time is 100 s, histograms similar 
to Fig. [7] are too noisy to be able to determine a rejecting 
criterion (i.e. a limiting highest acceptable value for and 
d\). Things become better only if one restricts to searching 
the brightest source, which is moreover supposed to be two 
magnitudes brighter than the second one. The other fainter 
star is not searched for by CLEAN. This case is still inter- 
esting because this configuration can occur at the galactic 
centre if a bright flaring source is observed when a faint star 
is also present in the field: only the flare's position would 
be of interest. Fig. [5] shows that the precision in the best 
direction io^\ m,K ^ 15 is here too of the order of 30 /xas. 
For all these 100 s results, the percentage of skipped runs 
is always zero, except for thk 16.75 where it remains less 
than 0.1 per cent. 

4.5 Three-star case 

Fig. [6] shows that when integrating the whole night, the 
precision in both directions is still of the order of the 
Schwarzschild ISCO for magnitudes higher than mm ~ 16. 

When the integration time is 100 s, and when searching 
only for the brightest source, the precision in the best direc- 
tion for magnitudes higher than mx = 15.5 is of the order 
of 40 fj,as, somewhat higher than the ISCO. The percentage 
of skipped runs is always less than 0.5 per cent. 

The main numerical results of this section are summa- 
rized in Table H 

All these results are very encouraging as far as the GR- 
related science cases of GRAVITY are concerned: the im- 
mediate vicinity of the black hole is within reach of the in- 
strument even if the source is quite complex. 



5 INFORMATION RETRIEVED WHEN A 
VARIABLE SOURCE IS IN THE FIELD OF 
VIEW 

In this section, the field of view contains a variable source 
with maximum magnitudes thr =13, 14 or 15. These values 
are realistic for the galactic centre flares: we have seen in 
Section [I] that the brightest flare observed to date has a 
maximum magnitude hik = 13.5. 

As for the brightness, it is modeled by the superimposi- 
tion of a Gaussian and a sine curve (see Fig. [sj> to take into 
account the two characteristic time-scales of the flares men- 
tioned in Section [l] It is assumed that the source light curve 
is measured with a 10 per cent accuracy by GRAVITY'S 

2 Here rriK is the magnitude integrated over the field of view, not 
the magnitude of the brightest star only. 
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Table 2. Imaging mode best precision as a function of integration time and sampling (ppm = pixels 
per mas), and percentage of skipped 'bad runs' (for 100 s integration). When more than one star 
are present and the integration time is 5 h, the stars are supposed to be more distant than the PSF 
characteristic dimension to avoid 'bad runs'. When the integration time is of 100 s, only the brightest 
source is searched. 



Nb of stars 5 h 100 s % skipped (for 100 s) 

1 13 /ias, 3 ppm 4 /ias, 10 ppm ^ 0.03 

2 25 /(as (not too close), 3 ppm 30 /ias (brightest), 3 ppm ^0.1 

3 30 /ias (not too close), 3 ppm 40 /ias (brightest), 3 ppm ^ 0.5 



adaptively corrected guiding cameras. This value of inten- 
sity will be used in the fitting procedure which will allow to 
obtain the source position. 

Two kinds of simulations are performed. Simulations of 
the first kind are done with a non-moving source, with vari- 
able magnitude. Simulations of the second kind are done 
with a source both variable in position and magnitude. Its 
trajectory is circular, at a distance of 30 /ias from the cen- 
tre, corresponding to the Schwarzschild ISCO of a black hole 
of 4.3 1O 6 M0 (in this simple model, the lensing effects dis- 
cussed e.g. in Hamaus et al.|2009 are neglected) . 

In all simulations, an observation is assumed every 100 s 
during the whole night. Each observation gives a set of visi- 
bilities and phases. The period of the night in which the ob- 
served intensity is higher than half the maximum observed 
intensity is selected. An image is then reconstructed for each 
set of visibilities and phases in this interval. Fifteen such 
100 s integration snapshots are thus obtained on both sides 
of the maximum of the flare. 

The astrometric positions measured on the recon- 
structed images are found by resorting to pure imaging mode 
(with a 10 pixels per mas sampling), or by fitting the com- 
plex visibilities with a single-source model using the imaging 
mode retrieved positions as first guesses. 

One possible way to check whether a motion can be 
detected is to calculate the dispersion of the positions thus 
obtained, in both directions: 



<5 X = 



8 V = 



\ 



14 



X) , 



(5) 



\ 



1 

14^ (j/k 



yf 



where x^ and j/ K are the retrieved positions of the flare in 
the x and y direction in each of the fifteen snapshots, and 
the barred quantities are the average position during the 
considered night. 

By making a Monte Carlo analysis based on this 
scheme, the histogram of the 2D dispersion can be com- 
puted. Fig. [9] represents the isocontours of these 2D his- 
tograms in the pure imaging mode case. Fig. [To] gives the 
same results when a fitting procedure is used. In both fig- 
ures, the blue contours are related to the simulations with a 
non-moving source, and the red contours are related to the 
simulations with a moving source. 

In practice, one night of observation will give one value 
of the dispersion in both directions (<5 x ,<5 y ). The interest of 



these histograms is that if this point falls into one of the 
two sets of contours (blue or red), and if the two sets do 
not overlap or very slightly overlap, then it is possible to 
conclude whether the flare is moving or not. 

Of course, the model used here is very simple, because 
no relativistic effects are taken into account to model the 
flare motion, and the orbit is supposed to be seen face on. 
The effect of the inclination is not taken into account: it 
would shift the histograms. Moreover, if the orbit is seen 
nearly edge-on, relativistic effects will make it look very 



different, as described extensively in Hamaus et al. (20091 



Consequently, Figs. [9] and |10| must be taken as a first step 
towards a more realistic modeling, which will be done in 
future work. 

It appears first on Fig. [9| and |10| that the histograms of 
non-moving sources (solid blue ones) are centered on a value 
of (S x , S y ) which is close to the value of the astrometric pre- 
cision (pure imaging mode precision, or 'imaging mode + fit' 
precision, respectively) corresponding to one point source of 
the same magnitude. Indeed the pure imaging mode astro- 
metric precisions for one star of magnitude ttik = 13, 14, 15 
are in /ias: (<r x = 7.7, o y = 4.8), (<r x = 13.6, cr y = 6.3) 
and (<r x = 28.6, (J y = 11.6). And the 'imaging mode + 
fit' corresponding precisions are: (<r x = 6.7, o y — 2.6), 
(er x = 12.8, cr y = 5.0) and (cr x = 27.4, a y = 10.6). More- 
over, it appears that the pure imaging mode results are as 
rich as the fitting ones: the histograms are very close to each 
other, which is not surprising as it has been shown in the 
preceding section that fitting does not allow to gain a lot in 
precision in the 100 s integration case. 

Secondly, it can be concluded that when the magnitude 
of the flare is brighter than 14, it is easy to distinguish be- 
tween the cases of a moving flare and a non-moving flare. 
For a flare fainter than ?tik = 15, it becomes more difficult 
to conclude because the histograms overlap. However, if we 
suppose that the flare is indeed due to a hotspot orbiting 
on the ISCO, it is possible to use the fact that the motion 
is not random and to fit the positions obtained throughout 
the night to the last stable orbit. The \ 2 given by this pro- 
cedure is given in Fig. [Tl] where the fit is performed for an 
mK = 15 flare moving on the ISCO, and for a non-moving 
flare: in the first case, the model is consistent with the data, 
whereas in the second case, the model is wrong. The question 
is to determine whether the precision on the flare positions 
is sufficient to distinguish the two cases. 

It appears that the two distributions are well separated: 
even for a flare of maximum magnitude tuk = 15, it is still 
possible to determine whether it is moving or not. 

However, if a flare of maximum magnitude m,K = 16 
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Flare light curve 
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Figure 8. Light curve of the variable source. Intensities are in arbitrary units and the abscissa is graduated in hours, UT. 



Dispersion histogram, no fit, mK=13 flare 



Dispersion histogram, no fit, itik=14 flare 



Dispersion histogram, no fit, mx=15 flare 




Dispersion in the x direction (pas) 
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Dispersion in the x direction (Mas) 
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Figure 9. Isocontours of the histograms of the 2D dispersion for a moving flare (dotted red) and a non-moving one (solid blue); the 
maximum magnitude of the flare is 13, 14, 15 for upper-left, upper-right and bottom pictures. The retrieved positions are found in the 
imaging mode. The axes are graduated in /tas. 



is considered, the histograms are impossible to distinguish, 
even if resorting to the ISCO fitting procedure. 

It is besides possible to estimate the typical size of the 
orbit provided the flare is bright enough. Indeed, in the sim- 
ulation, the flare makes one orbit in 20 min, i.e. it takes 12 
snapshots of 100 s each to cover a full orbit. The angle cov- 
ered between two successive snapshots is of 6 — y|. Then 
the theoretical dispersion (the same in both directions) is, 
with .Risco being the ISCO radius, and keeping in mind that 
15 snapshots are simulated in total: 



\ 



g X>L„™» 2 (") 

fc=0 



(6) 



23 /ias 



which is quite close to the centre of the red dotted histogram 
in Fig.[9|and |10| for tuk =13 and 14. Thus, the characteristic 
size of the orbit is measurable if the flare is as bright at its 
maximum as dik = 14. 



6 CONCLUSION 

We have analysed the astrometric information that an in- 
strument such as GRAVITY is able to provide when used 
in 'imaging' mode. The accuracy of the positions measured 



for a field of fixed stars has been investigated. A variable 
source orbiting on the ISCO of Sgr A* has then been con- 
sidered. 

The first study leads to the conclusion that the astro- 
metric precision in pure imaging mode is of the order of, 
or better than the Schwarzschild radius of the black hole 
when a single bright enough star is in the field. The limit- 
ing brightness of the source is a function of the integration 
time: rriK < 16 is needed for the whole night of integration to 
reach this accuracy, whereas rriK < 15 is needed when the 
integration time is 100 seconds. These results lead to the 
conclusion that the imaging mode already provides valuable 
astrometric performance. Besides, when used to bootstrap a 
fitting procedure, the precision reaches the theoretical limit 
of the instrument. 

The second study makes it clear that it is possible to 
distinguish between a moving and a non-moving flare on 
the ISCO, even for flares with a maximum magnitude of 
only rriK = 15. The positional information can be retrieved 
in the imaging mode (i.e., without any fit, so without any 
model). Moreover, the characteristic size of the orbit is at 
hand provided the flare is as bright as dik = 14. 

These results are comforting and will be developed in 
future work. The model of the source is quite crude for the 
time being and must be refined: for this purpose, a general 
relativistic ray-tracing code is needed, as well as a more re- 
alistic model for the flare. Still, all these results are very 
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Dispersion histogram, with fit, mi<=13 flare 



Dispersion histogram, with fit, mK=14 flare 



Dispersion histogram, with fit, mK=15 flare 
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Figure 10. Same as Fig. [9] the retrieved positions being found by resorting to fitting, the imaging mode results being used as a first 
guess. 



X 2 histograms for an ISCO fit with an m K =15 flare 




Figure 11. Histograms of reduced x 2 f° r a non-moving flare (blue) and a flare orbiting around the ISCO (red). This \ 2 is the root mean 
square of the difference between the retrieved positions and a model assuming a circular motion with radius 30 ^tas. 



encouraging for the GR-related science cases of the instru- 
ment: the immediate vicinity of the black hole horizon is 
within reach. Having access to the motion of the flare would 
also allow constraining the black hole spin. Indeed, the tra- 
jectory of the hotspot is a function of this parameter. How- 
ever, this constraint will only be possible provided an ade- 
quate hotspot model is available. 
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Hence the shot noise: CT s hot = y (m) . It is assumed that 
the number of read out noise electrons is equal to £>r ON = 36. 
Given the rate of dark current electrons iVdark = 100 s , 
the variance of the detection noise per integration time r is: 

Cdetec = "'shot + "RON + ^dark T ■ (A5) 

Taking this noise into account, and assuming a gaus- 
sian distribution for the various noise contributions, it is 
possible to compute a realization of the detection noise 
riA, n,B, nc, no corrupting the four A,B,C,D signals. The sig- 
nal to noise ratio is lower for fainter sources: it is magnitude 
dependent. 



APPENDIX A: MODELLING THE 
INTERFEROMETRIC OBSERVABLES 

GRAVITY will use the four VLT UTs to observe Sgr A*. 
Each baseline will combine the beams of intensities I p and I q 
delivered by two telescopes labeled p and q. Let the intrinsic 
visibility modulus and phase of the observed object be Vobj 
and <f>obj- Let S(f> be the phase shift between the two chan- 
nels from telescope p and q. The short-exposure noiseless 
combined intensity is: 



1(6$ 



!> + 0obj+^piston, S hort)(Al) 



where <^ p iston, short is the short-exposure piston phase. 

It is sufficient to get four samples of this function in 
order to retrieve the visibility modulus and the phase of the 
object. Indeed, if the four following quantities are computed: 



A 
B 

C = /(- 

D = 



1(0), 

iC-) 



(A2) 



then it is easy to show that the complex visibility of the 
object is given by: 



Vobje 



i^ = M-C)-i{B-D) 
A+B+C+D ' 



(A3) 



In order to be able to simulate realistic quantities, we 
now take into account the various sources of noises affecting 
the different quantities that have been introduced. 



Al Detection noise on the four samples A, B, C 
and D 

Let iVph be the number of photons arriving from each of 
the JVtai telescopes of the interferometer. Each set of iVph 
photons will be dispatched to the JVtei — 1 other telescopes. 
Then N p _^ 1 photons will be present in the baseline between 
two telescopes. Let T be the transmission of the instrument 
multiplied by the quantum yield, estimated to be 0.009. The 
mean number of photo-electrons (m) per sample per baseline 



(m) 



2jV P hT 
4(A tol - 1) ' 



(A4) 



A2 Strehl ratio fluctuation 

It must also be taken into account that the AO-corrected 
intensities I p and I q are noisy because of Strehl ratio fluc- 
tuations. Let n v and n q be realizations of these noises. 



A3 Visibility modulus and phase noise 

One determination of the visibility modulus and phase us- 
ing equation ( |A3[ ) is made after having integrated during 
100 s. One value of the complex visibility obtained in this 
way can thus be seen as an average of short exposure com- 
plex visibilities which would be obtained after an integration 
time smaller than the turbulence evolution time (of the or- 
der of the ms). Each short exposure complex visibility can 
be written: 



V s 



Vobje 



i (0obj+0piston, short) 



(A6) 



where </> p iston, short is the piston phase, which varies at a rate 
equal to the evolution time of the turbulence. 

Averaging these short exposure quantities, the 100 s 
integrated complex visibility becomes: 



(A7) 



where C is a contrast parameter between and 1, and </>pi s ton 
is the result of the averaging of the different short exposure 
piston phases. This is the complex visibility obtained at the 
science beam. 

In order to get rid of the piston phase, GRAVITY will 
observe both the science object and a reference star of zero 
intrinsic phase within the same isopistonic patch. The in- 
strument fringe tracker stabilizes the fringes of the reference 
and applies the same correction to the object: the science 
fringe pattern is thus stabilized. However, the two sources 
being separated on the sky, this procedure is affected by an 
error depending on the angular separation 6 of the stars: a 
piston residual thus remains. This anisopistonism error on 
the location of the source on the sky can be computed for 



a typical Paranal atmosphere using Delplancke (2008) and 



Shao & Colavita ( 1992 1 



^piston res [arcsec] = 370 



Lo[m]3 0[rad] 



(A8) 



where B is the baseline, Lq is the outer scale of the turbu- 
lence and T the integration time. 

This allows to compute the residual piston phase due 
to anisopistonism: 
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-"piston res 
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piston res 5 



(A9) 



where Q p i s ton res is a realization of the anisopistonism error 
on the location of the source on the sky, and A is the obser- 
vation wavelength. 

Moreover, the stabilization of the science fringe pattern 
on the reference position is done within a certain precision. 
A residual phase due to the limited performance of the fringe 
tracking must also be taken into account. It will be denoted 

^tracking res ■ 

Finally, the origin of angular positions on sky is by con- 
struction the location of the reference star. The phase of the 
object is thus shifted by me t = ^AS • B, where B is the 
vector joining the two telescopes (baseline vector) and AS 
is the difference between the position vector of the scientific 
and reference targets on the sky. This shift is determined by 
the instrument's metrology system. It is affected by an error: 
^met res = 0.07 rad for an integration time of 30 ms, which 
is one of the design parameters for GRAVITY (Ra bien et al.| 
2008b. 



Putting all this together, it becomes possible to express 
the observed values (denoted with a hat) of the visibility 
modulus and phase on the science beam: 



C) 2 + (B~D)' 2 



(A + B + C + D) 2 
and the phase is: 



4°d ctec 



= at an 



D-B 
A-C 



(A14) 



(A15) 



These values of visibility and phase are delivered by 
a program that simulates the instrument GRAVITY. The 
noises affecting these quantities are thus at hand, and will 
be used for the computation of the image from the interfer- 
ometric observables (see Section [jjj. 

This paper has been typeset from a TgX/ P/TjrjX file prepared 
by the author. 



V BC = CV ohi +n Vs<: , (A10) 

0sc — 0obj ^piston res ~t~ ^tracking res ~t~ 0mct res ~\~ ^tp sc , 

where ny sc and ri0 ac are realizations of the noises affecting 
the measures of the visibility and the phase. 



A4 Effective computation in the simulations 



The long-exposure equivalent of equation (All which takes 
all noise sources into account is: 



I{5<t>) = I p + n p + I q + n q 



(All) 



+ 1\J {I P + n p )(I q + n q ) V sc < 

~\~ ^detec; 



where ridctcc is the realization of the detection noise made 
explicit in Section |Al| 

This quantity can be computed for each of the four sam- 
ples: 

A = /((», (A12) 

6 . /(I), 

C = />), 

b = /(^), 

Finally the simulation of the measured complex visibil- 
ity is computed according to: 

A+B+C+D 

In order to work with de-biased quantities, the squared 
visibility is computed according to: 



